0. Origin of the Neural Networks : the perceptron

A neurology primer :

A neuron is a cell which can be electrically excited. It processes and transmits information in the brain. Typically, a neuron is made of a cell body, dendrites, and an axon. Networks of neurons are interconnected in the brain by membranes called synapses. Synaptic signals from other neurons are received at the dendrites while signals to other neurons are transmitted by the axon.

Perceptrons were first introduced in the 1950s and 1960s as an artificial neuron.

It is easy to see the analogy between the two : in an artificial neuron (perceptron), several inputs are combined (mimicking the role of the dendrites) to produce an output (mimicking the role of the axon).

A simple way to combine the inputs into a single output is to compute a weighted sum. The value of this sum with respect to a given threshold gives the output.

\begin{eqnarray} \mbox{output} & = & \left\{ \begin{array}{ll} 0 & \mbox{if } \sum_j w_j x_j \leq \mbox{ threshold} \\ 1 & \mbox{if } \sum_j w_j x_j > \mbox{ threshold} \end{array} \right. \end{eqnarray}

We an actually replace the threshold by introducing what we call a bias term :

$$\begin{eqnarray} \mbox{output} = \left\{ \begin{array}{ll} 0 & \mbox{if } w\cdot x + b \leq 0 \\ 1 & \mbox{if } w\cdot x + b > 0 \end{array} \right. \end{eqnarray}$$

This will prove handy in the following.

It is quite easy to visualise the perceptron as a way of weighing information. Assuming we want to answer the question : will it rain tomorrow ?. We can gather some data about the weather today :

• Temperature today
• Wind direction
• Weather in a 100 km radius

The perceptron combines all this information with varying importance (the weights) and eventually fires a yes or no answer.

The weights can be determined using existing data and minimising the prediction error on this data with gradient descent.

Why stop at a single neuron ? After all, the brain is made up of about 100 billion neurons! We can combine several perceptrons together (see fig) in the hope of creating greater levels of abstractions which will allow for more complex decisions.

In the following, we will now see how the weights and biases in a complicated network can be tuned to adapt automatically to any data.

1. From the perceptron to the sigmoid neuron

Recall that perceptrons take decisions on hard thresholding. This makes it difficult to tune our network. Let's give an example. Sat we want to improve the network's performance on misclassified data points. We can make a small change on one of the weight and see if classification improves. However, with hard thresholding, small changes may completely alter the behaviour of the network as the output may brutally fall on the other side of the threshold. This problem is alleviated by the use of a continuous activation function, like the sigmoid (soft thresholding).

With the sigmoid, the output of the neural network becomes :

$$\begin{eqnarray} \frac{1}{1+\exp(-\sum_j w_j x_j-b)}. \end{eqnarray}$$

The differences between hard and soft thresholding are illutrated below

In fact, any function with a shape similar to the sigmoid (a smooth transition from 0 to 1 between) would suit our need of soft thresholding. It just so happens that the sigmoid has properties which make it particularly attractive.

2. Architecture of the neural network

Now that we have the base unit of our network figured out, how do we build the full network ? There are arbitrarily many complicated ways of building a network (think feedback loops, arbitrary connections from one layer to another ...). We are going to focus on a subclass of networks, which will prove handy to parameterize : feedforward fully connected neural networks.

Such networks are organised in layers, each unit of each layer being connected to each unit of the following layer.

3. Training the neural network : backpropagation

From now on, we'll work with the MNIST data set to provide an illustration to our discussion. This dataset is organised in the following way: It is split into two parts. There is a training part (60 000 images, represented as 28*28 arrays with the grayscale value of the relevant pixel) and a test part (10 000 similar images). Each image corresponds to a single handwritten digit.

Architecture :

We will use a simple design : 3 layers (inputs + one hidden layer + one output layer).

One such network is shown below, with 15 hidden layer cells

Notations :

x is the input, a 28*28 = 784 dimensional vector y is the output, a 10-dimensional vector such that (1,0,0,0,0,0,0,0,0,0) represents 0, (0,1,0,0,0,0,0,0,0,0) represents 1 and so on. a is the 10-dimensional output of the neural network

We will quantify the success of our model with the well known squared loss error :

$$\begin{eqnarray} C(w,b) \equiv \frac{1}{2n} \sum_x \| y(x) - a\|^2. \end{eqnarray}$$

The reason behind this choice is that the square loss function is smooth in the parameters of the networks (i.e. its weights and biases) which makes it easier to understand how small changes in the weights and biases may affect the classification performance.

Gradient descent :

In order to minimise the error $C$ identified above, we are going to use gradient descent.
In fact, we will make use of a slightly modified version of gradient descent to alleviate computing constraints. Indeed, vanilla gradient descent updates the gradient of $C$ by computing the derivatives with respect to all data points at each update. If we have a large training set, this may prove untractable.
That's why we'll implement stochastic batch gradient descent :

• We will randomly select some data points: a mini batch
• And we will compute the gradient update with respect to these points only
• We continue sampling the training set randomly (without replacement) until it is exhausted. This completes an epoch of training.
• We iterate over several epochs.

The size of the batch can vary, with a small size meaning fast computations at the cost of a loss of accuracy in the computation of the gradient.

Let us recapitulate : our goal is to tailor the weights $w$ and the biases $b$ of the neural network so that we minimise the cost $C$.
Gradient descent helps us achieve this goal by providing a way to progressively update all $w$ and $b$ to their ideal values.

Mathematically :

$$ \begin{eqnarray} w & \rightarrow & w' = w-\frac{\eta}{m} \sum_j \frac{\partial C_{X_j}}{\partial w} \\ b & \rightarrow & b' = b-\frac{\eta}{m} \sum_j \frac{\partial C_{X_j}}{\partial b}, \end{eqnarray}$$

where the sum is on a mini-batch.

Backpropagation

Let's move on to the more difficult bit.
The backpropagation algorithm provides us with an efficient way of computing the gradient of $C$. Let's see how.

Notation

$b^l_j$ is the bias of the jth neuron in the lth layer.
$a^l_j$ is the activation of the jth neuron in the lth layer.
$w^l_{jk}$ is the weight for the connection from the kth neuron in the (l−1)th layer to the jth neuron in the lth layer.

Let us write down a first equation to write the output of a layer :

\begin{eqnarray} a^{l}_j = \sigma\left( \sum_k w^{l}_{jk} a^{l-1}_k + b^l_j \right), \end{eqnarray}

This equation expresses the architecture of our network: the sum is on all the neurons of the previous layer, thereby making our network fully connected.
We will use the following shorthand notation as well :

$$z^l_j = \sum_k w^l_{jk} a^{l-1}_k+b^l_j$$

This immediately gives :

$$a^l = \sigma(z^l)$$

Computing the derivatives of $C$ :

Now that the notations are in place, recall that we want to compute the derivatives of $C$: $\hspace{2mm}\partial C / \partial b^l_j$ and $\partial C / \partial w^l_{jk}$.

a) Computing $\hspace{2mm}\partial C / \partial b^l_j$ :

Using the chain rule, we write :

\begin{equation} \frac{\partial C}{\partial b^l_j} = \sum_k \frac{\partial C}{\partial z^l_k}\frac{\partial z^l_k}{\partial b^l_j } \end{equation}

which reduces to :

\begin{equation} \frac{\partial C}{\partial b^l_j} = \frac{\partial C}{\partial z_l^j} \equiv \delta^l_j \end{equation}

where we have introduced the intermediate quantity $\delta^l_j \equiv \frac{\partial C}{\partial z_l^j}$.

b) Computing $\delta^l_j $ :

Let's compute it in two times.

For the final layer $L$, using the chain rule,

\begin{eqnarray} \frac{\partial C}{\partial z^L_j} &= \sum_i \frac{\partial C}{\partial a^L_i}\frac{\partial a^L_i}{\partial z^L_{j} }\\ &= \frac{\partial C}{\partial a^L_j}\frac{\partial a^L_j}{\partial z^L_{j}} \\ & = \frac{\partial C}{\partial a^L_j} \sigma^{\prime}(z^L_j) \end{eqnarray}

For any other layer $l$, using the chain rule:

\begin{eqnarray} \delta^l_j & = & \frac{\partial C}{\partial z^l_j}\\ & = & \sum_k \frac{\partial C}{\partial z^{l+1}_k} \frac{\partial z^{l+1}_k}{\partial z^l_j}\\ & = & \sum_k \frac{\partial z^{l+1}_k}{\partial z^l_j} \delta^{l+1}_k \end{eqnarray}

note that :

\begin{eqnarray} z^{l+1}_k & = & \sum_j w^{l+1}_{kj} a^l_j +b^{l+1}_k\\ & = & \sum_j w^{l+1}_{kj} \sigma(z^l_j) +b^{l+1}_k \end{eqnarray}

hence :

\begin{eqnarray} \frac{\partial z^{l+1}_k}{\partial z^l_j} = w^{l+1}_{kj} \sigma'(z^l_j). \end{eqnarray}

thus :

\begin{eqnarray} \delta^l_j & = \sum_k w^{l+1}_{kj} \delta^{l+1}_k \sigma'(z^l_j) \\ \delta^l & = ((w^{l+1})^T \delta^{l+1}) \odot \sigma'(z^l) \nonumber \end{eqnarray}

The last line expresses the result using the Hadamard product.

c) Computing $\partial C / \partial w^l_{jk}$ :

All that's left is to compute $\partial C / \partial w^l_{jk}$. Once again, we make use of the chain rule :

\begin{eqnarray} \frac{\partial C}{\partial w^l_{jk}} & = \sum_i \frac{\partial C}{\partial z^l_i}\frac{\partial z^l_i}{\partial w^l_{jk} }\\ & = a^{l-1}_k \delta^l_j \end{eqnarray}

The backpropagation algorithm

The formulas we established guide us towards the backpropagation algorithm :

• Initialise by setting $a^1$ as the input data
• Feedforward : propagate by computing $z^l$ and $a^l$ for each leayer
• Compute $\delta^L$.
• Backpropagate : compute $\delta^l$ for all remaining layers
• Use $\delta$ to compute $\hspace{2mm}\partial C / \partial b^l_j$ and $\partial C / \partial w^l_{jk}$.

4. Python implementation

The code below is an implementation of neural networks trained by backpropagation in pure python.

As advertised above, we use the sigmoid as our activation function and the mean squared error as our cost function. The other important parameters are :

• sizes : specifies for each layer the number of unit. Here, we have 3 layers. The input layer has 784 units, the hidden layer has 30 units and the output layer has 10.
• batch_size : The number of samples in the mini batch
• epochs : The number of training epochs
• lambda : The regularisation parameter (this is more advanced, stick to 0)
• learning_rate : The learning rate in the gradient descent algorithm. This quantifies how quickly we follow the direction of increasing gradient

You can now play around with the parameters and see whether you can improve the classification accuracy ! (Change number of layers, number of units, learning_rate ...)